home *** CD-ROM | disk | FTP | other *** search
/ io Programmo 60 / IOPROG_60.ISO / soft / c++ / gsl-1.1.1-setup.exe / {app} / src / rng / fishman2x.c < prev    next >
Encoding:
C/C++ Source or Header  |  2001-12-02  |  2.5 KB  |  112 lines
  1. /* rng/fishman2x.c
  2.  * 
  3.  * This program is free software; you can redistribute it and/or modify
  4.  * it under the terms of the GNU General Public License as published by
  5.  * the Free Software Foundation; either version 2 of the License, or (at
  6.  * your option) any later version.
  7.  * 
  8.  * This program is distributed in the hope that it will be useful, but
  9.  * WITHOUT ANY WARRANTY; without even the implied warranty of
  10.  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11.  * General Public License for more details.
  12.  * 
  13.  * You should have received a copy of the GNU General Public License
  14.  * along with this program; if not, write to the Free Software
  15.  * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
  16.  */
  17.  
  18. /*
  19.  * This generator is taken from
  20.  *
  21.  * Donald E. Knuth
  22.  * The Art of Computer Programming
  23.  * Volume 2
  24.  * Third Edition
  25.  * Addison-Wesley
  26.  * Page 108
  27.  *
  28.  * It is called "Fishman - L'Ecuyer"
  29.  *
  30.  * This implementation copyright (C) 2001 Carlo Perassi.
  31.  */
  32.  
  33. #include <config.h>
  34. #include <stdlib.h>
  35. #include <gsl/gsl_rng.h>
  36.  
  37. /* Fishman */
  38. #define AAA_F 48271UL
  39. #define MMM_F 0x7fffffffUL    /* 2 ^ 31 - 1 */
  40.  
  41. /* L'Ecuyer */
  42. #define AAA_L 40692UL
  43. #define MMM_L 0x7fffff07UL    /* 2 ^ 31 - 249 */
  44. #define QQQ_L 52774UL
  45. #define RRR_L 3791UL
  46.  
  47. static inline unsigned long int ran_get (void *vstate);
  48. static double ran_get_double (void *vstate);
  49. static void ran_set (void *state, unsigned long int s);
  50.  
  51. typedef struct
  52. {
  53.   unsigned long int x;
  54.   unsigned long int y;
  55.   unsigned long int z;
  56. }
  57. ran_state_t;
  58.  
  59. static inline unsigned long int
  60. ran_get (void *vstate)
  61. {
  62.   ran_state_t *state = (ran_state_t *) vstate;
  63.  
  64.   long int y = state->y;
  65.   long int r = RRR_L * (y / QQQ_L);
  66.  
  67.   y = AAA_L * (y % QQQ_L) - r;
  68.   if (y < 0)
  69.     y += MMM_L;
  70.  
  71.   state->y = y;
  72.   state->x = (AAA_F * state->x) & MMM_F;
  73.   state->z = (state->x - state->y) & MMM_F;
  74.  
  75.   return state->z;
  76. }
  77.  
  78. static double
  79. ran_get_double (void *vstate)
  80. {
  81.   ran_state_t *state = (ran_state_t *) vstate;
  82.  
  83.   return ran_get (state) / 2147483648.0;
  84. }
  85.  
  86. static void
  87. ran_set (void *vstate, unsigned long int s)
  88. {
  89.   ran_state_t *state = (ran_state_t *) vstate;
  90.  
  91.   if (s == 0)
  92.     s = 1;            /* default seed is 1 */
  93.  
  94.   state->x = s & MMM_F;
  95.   state->y = s & MMM_L;
  96.   state->z = s & MMM_F;
  97.  
  98.   return;
  99. }
  100.  
  101. static const gsl_rng_type ran_type = {
  102.   "fishman2x",            /* name */
  103.   MMM_F,            /* RAND_MAX */
  104.   0,                /* RAND_MIN */
  105.   sizeof (ran_state_t),
  106.   &ran_set,
  107.   &ran_get,
  108.   &ran_get_double
  109. };
  110.  
  111. const gsl_rng_type *gsl_rng_fishman2x = &ran_type;
  112.